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1 Introduction 

To construct orthonormal compactly supported wavelets we use as the starting point the 
dilation equation 

(j)(x) = v / 2^c fc 0(2x - k), i6l, c fe GM, (1) 

for the scaling function or father function 0. A solution exists in the sense of distributions 
when ^ k Ck = \/2 and it has compact support if only finitely many of the filter coefficients 
Cfc are non-zero (in fact, if Ck ^ only for k G {0, . . . , n} then the support of <fi is contained 
in [0, n]). The solution of the dilation equation is given by 

oo 

fe)=n m °( 2 ^)' ( 2 ) 
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where the trigonometric polynomial mo(£) is the scaling filter of <p, 

m ^) = ^5> e ~^ ( 3 ) 

The condition when the solution to (0) yields wavelets can be expressed in terms of 
the scaling filter m (£). If m satisfies 

|m (0| 2 + l^o(vr + 0| 2 = l, ^M, (4) 

m (0) = 1, and a technical condition called the Cohen criterion holds, then there exists 
a scaling function and wavelet pair Jl|, 0, Section 6.3]. We will not need exact details 
about Cohen's condition, we only use that if mo does not vanish on [— tt/2, tt/2], then it 
satisfies the criterion. 

If (f> is such that it gives rise to a wavelet, the wavelet ip is obtained by 

i>(x) = Y,(-V k Ci-k<P(2x-k). 

kez 

Note that then and i[) have the same regularity properties. 

The original construction of orthonormal compactly supported wavelets by 
I. Daubechies || considers a trigonometric polynomial mo with at most 2N non-zero 
coefficients c , . . . , c 2 n-i- The polynomial m is required to satisfy (Q) and to have a 
zero at tt of maximal order N (the order of the zero is also the number of vanishing 
moments the corresponding wavelet has). The resulting equations can be solved and the 
constructed mo has no roots on (— tt, tt), hence it satisfies the Cohen criterion. (Q is a 
comprehensive reference for these results.) 

To construct smoother wavelets our approach is to reduce the order of the zero m has 
at ir and instead to introduce zeros on the interval (tt/2, tt). Keeping N and the number 
of zeros on (tt/2,tt) fixed we use numerical optimization to choose the locations of the 
roots so that the resulting wavelets are as smooth as possible. We are able to construct 
wavelets that are more regular than those introduced in ||, M, and [ |I~5|| . 

This paper is organized as follows: in Section ^| we state our results and compare them 
to the above mentioned papers. The methods used are described in Section and in 
Section f| the optimization is discussed in detail. We comment in Section |5| on some other 
approaches that were tried. Finally, filter coefficients for selected most regular wavelets 
are listed in section |6]. 

2 Results 

In the numerical optimization for smoothness we studied wavelets with up to 40 filter 
coefficients and whose scaling filter has up to four roots on the interval (tt/2, it). The 
Sobolev regularity exponents of the most regular wavelets found are shown in Table |l|. 
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1.00 
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4.00 


2.00 
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3.64 


4.42 


4.73 


4.76 


3.00 


24 


3.87 


4.78 


5.14 


5.22 


4.00 
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4.11 


5.14 


5.55 


5.67 


5.00 
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4.34 


5.50 


5.95 


6.10 


5.84 


30 


4.57 


5.85 


6.34 


6.51 


6.38 


32 


4.79 


6.19 


6.71 


6.89 


6.88 


34 


5.02 


6.52 


7.09 


7.25 


7.30 


36 


5.24 


6.83 


7.45 


7.62 


7.69 


38 


5.47 


7.15 


7.81 


7.99 


8.08 


40 


5.69 


7.46 


8.17 


8.37 


8.51 



Table 1: The best Sobolev regularity exponents found. The length of the filter is 2N and 
the filter mo has n z roots on the interval (ir/2,ir). For comparison the column n z = 
gives the Sobolev exponents of the Daubechies wavelets ||, the other columns are new. 

(By the Sobolev regularity exponent s of <p we mean s = sup{s: <p G H s }, where H s is 
the Sobolev space consisting of all / G L 2 (R) such that (1 + |£| 2 ) s/2 /(0 G L 2 (R).) The 
results are also summarized in Figure 

The column n z = gives the Sobolev exponents for the Daubechies wavelets con- 
structed in [[| . For N > 4 we have constructed wavelets which are more regular. For 
large values of iV there are four families of wavelets corresponding to different numbers 
of roots on (vr/2, n), and the regularity increases with the number of roots. 

To relate these results to Holder regularity, it can be shown that the Holder regularity 
exponent a Q = sup{a: G C a } satisfies a G [s — 1/2, Sq] |[3|, Corollary 9.9]. Here C a 
consists of functions that are continuously different iable up to order equal to the integer 
part of a, and the highest order continuous derivative belongs to the (standard) Holder 
class with exponent equal to the fractional part of a. 

Note that in each column n z — 1, . . . 4, the first few entries are smaller than in the 
previous columns on the same row. This is not surprising, since it can be shown that in 
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Figure 1: The Sobolev regularity exponents Sq from Table |l|. For clarity only those are 
shown that are better than the ones with the same length filter but with fewer roots on 

(7r/2,7r). 



order to have ip G H n , the wavelet ip must have at least n + 1 vanishing moments 13, 



Proposition 8.3]. The wavelets corresponding to the entries in the table have N — 2n z 
vanishing moments. 

We show below (see Section |4.2| ) that the entries in the column n z = 1 are best possible 
in the sense that the corresponding wavelets are smoothest possible among those with 
2iV > 8 coefficients and at least N — 2 vanishing moments. 

For higher values of n z the exponents listed in Table [l] are best possible among those 
wavelets with 2N coefficients and whose scaling filter m has the indicated number of 
zeros on (7r/2,7r). Some of these were subjected to extensive tests by choosing many 
different initial approximations for the optimizer. Note also how the results, except for 
the first few entries, lie on almost straight lines in Figure |]. See also Figures |3| and [7] 
which show typical examples of the regularity exponent as a function of the location of 
the roots. 

However, it is not completely clear if the results for n z > 2 are best possible if we 
consider wavelets with the same number of vanishing moments as above but do not insist 
on m having roots on (7r/2, 7r). We come back to this at the end of Section |4~3 . 

We remark that the case 2iV = 10, n z = 1 appears already in (TO . 



More regular wavelets than the original Daubechies family have been constructed by 
I. Daubechies 0, H. Volkmer [jT5 |, P. G. Lemarie-Rieusset and E. Zahrouni ||. We 



discuss the wavelets from H in Section 4.1 
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H. Volkmer [15] computes only the asymptotic regularity ratio \miN^ 00 a/2N and 
constructs wavelets for which the ratio is close to 0.19. Note that for many entries in 
Table [L] the ratio a/2N is over 0.21. 

P. G. Lemarie-Rieusset and E. Zahrouni || discuss what they call the Matzinger 
wavelets: the scaling filter has a multiple root at 2ir/3 and no other roots on (tt/2,tt). 
Our results are better, e.g., the smoothest one they construct has 34 coefficients and 
Sobolev exponent 6.13. Table |T| shows that already with n z = 1 our wavelet with the 
same number of coefficients is more regular. In fact, for each value of N the n z = 1 
wavelet listed in Table |l| is better than the smoothest wavelet discussed in f8j. 



3 Methods 

We begin by describing the theorem that allows us to compute exactly the Sobolev 
regularity exponent of solutions to ([]]). Define another trigonometric polynomial, r(£), 
by factoring out all zeros of |mo(£)| 2 at £ = 7r, that is, let 

K(f)l 2 =(^) M r (fl, (5) 

where M G N is chosen so that r{ji) ^ 0. The function r is used to define an operator 
on C[0,n] by 

T r u(C) = r(£/2)u(£/2) + r ( 7r _ ^ 2 )u(ti - £/2). (6) 

The Sobolev regularity exponent of a solution to the dilation equation ([I]) is given by the 
following theorem, proved independently by L. Villemoes and T. Eirola: 

Theorem 3.1 (|]7|, |T3]]) Suppose <fi is a solution to ([I]), M , r , and T r are as above, and 
r satisfies the Cohen criterion. Then <p G H s if and only if s < sq = M — log 4 p{T r ), 
where p(T r ) is the spectral radius ofT r on C[0,7r]. 

If r is a polynomial of cosines of degree at most d, i.e., r(£) = J2k=o ^ cos(k^), then 
T r maps the space spanned by {cos(k^)}f =0 into itself and the spectral radius of T r can 
be calculated from the eigenvalues of a finite matrix 0, O . 

Suppose we want to look at a filter of length 2N (i.e., a filter mo with at most 2N 
non-zero coefficients), such that mo has a zero of order M at it and n z zeros on the 
interval (7r/2,7r). Note that M is also the number of vanishing moments the associated 
wavelet has (if it exists). 

Define a new sequence of coefficients by 



2N-1 
k=0 



(7) 
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The orthonormality condition (^) is equivalent to 

a = -, a 2fc = 0, for k = 1, . . . , N - 1. 

We first use the zero of m at 7r to find (linear) relations among the {a^}: we express a^, 
fc = 1, 3, . . . , 2M - 1, in terms of a k , k = 2M + 1, 2M + 3, . . . , 2N - 1. Note that 
there are N — M = 2n z independent parameters. For convenience we switch to a new set 
of parameters v = (v\, . . . , v Hz ) that are simply obtained from ak, k = 2M + 1, 2M + 3, 
. . . , 2N — 1 by scaling the latter so that v = (1, . . . , 1) corresponds to the Daubechies 
wavelets (i.e., all roots of m are at it). The function r is then obtained by factoring 
((1 + cos£)/2) from |mo(£)| 2 . (The computations are done quickly with a symbolic 
algebra system such as Mathematica by using the Taylor polynomial of mo at it.) 

When r has low degree it is easy to compute the matrix of the operator T r on the 
span of {cos(/c£)}^!f r algebraically — an example is done in Section fO| . When the degree 



is high it is more convenient to compute the matrix numerically by using the fast Fourier 
transform. Finally, sq is computed from the largest eigenvalue of this matrix. 

Note that fixing the value of v does not uniquely determine the wavelet, since it only 
determines the absolute value of m (£)- For our purposes this is enough, since all those 



wavelets have the same regularity by Theorem 3.1 



4 Optimization 

4.1 One degree of freedom 

First consider wavelets with four filter coefficients, i.e., N = 2, and one vanishing moment. 
Then 

l m o(0| 2 = 2 + fll cos £ + a 3 cos 3£ 
and m(7r) = gives a% = 1/2 — a 3 . To compute r we get 

K(0| 2 = ^(1 + c °s 0(8a 3 cos 2 e-8a 3 cos£ + 1) = i(l + cosOr(0- 

The Daubechies wavelet has two vanishing moments, hence r(ir) = 0, which implies 
a 3 = —1/16. Introduce the new parameter v by a 3 = —v/16. Then 

1) 1) 

r (0 = ~2 cos2 £+ 2 cos £+ !> 

and the matrix of the operator T r on the span of {1, cos£, cos 2^} is easily computed to 
be 




8 - 2v -2v \ 
2v 2v 
— v 8 — 2v —v J 
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2N Hdlder Sobolev 



4 0.55 1.00 

6 1.09 1.42 



Table 2: Hdlder and Sobolev regularity exponents for the wavelets constructed in [[J. 

2N HSlder Sobolev 



4 0.59 0.87 

6 1.40 1.31 



Table 3: Holder and Sobolev regularity exponents for the Holder smoothest wavelets 
constructed in |J. (The Holder exponents are lower bounds for the actual value.) 

The eigenvalues of this matrix are 

Ai = 1 + \/l — v , A 2 = 1 — y/l — v, and A 3 = — — . 

Recall that r(£) has to be non-negative. In particular r(7r) > gives v < 1. Then Ai > 1 
and p(T r ) > 1 and so So < 1 for all admissible values of v. The Daubechies wavelet 
corresponds to v — 1 and So = 1, for all other wavelets in this family we have v < 1, and 
therefore s < 1. 

Numerical results show this is also the case for wavelets with larger values of N and 
with one degree of freedom. Thus the original Daubechies wavelets are the smoothest 
possible among wavelets with 2N coefficients and at least N — 1 vanishing moments. 

I. Daubechies || has also constructed wavelets that have higher Holder regularity 
exponents than the ones originally constructed in ||. In || the cases 2N = 4 and 
2N = 6 were studied. Table § gives both Holder and Sobolev exponents for two wavelets 
from the original construction in . (The Sobolev exponents are as in Table [I], the Holder 
exponents are from 0.) Table [3] gives the exponents for the two wavelets constructed 
in ||. The Sobolev exponents for the new wavelets are easily computed from Theorem 3.1 , 
e.g., for 2N = 4 we get s = log 4 (10/3) w 0.87. 

Comparing the numbers in Tables ||] and |3| an interesting phenomenon appears: both 
for 2iV = 4 and for 2iV = 6 we have an example of two functions, one of which is slightly 
more regular when regularity is measured in the Holder sense, but the opposite is true if 
Sobolev regularity is used. 

4.2 Two degrees of freedom 

Optimization was done on the v = (v\, v 2 ) parameters. Since some values of v yield an r v 
that takes on also negative values, this is now a constrained optimization problem with 
the nonlinear constraint mino<,t< 7r r v (£) > 0. 
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Figure 2: The Sobolev regularity exponent so as a function of the parameters v = (vi,V2) 
for the case 2N = 10. (When v is outside the admissible region minr v > we plot sq as 
zero.) 

Figure ^| shows a plot of the Sobolev exponent as a function of the parameters v±, Vi- 
A typical example of the exponent on the boundary of the feasible region is given in 
Figure |3], where the exponent is plotted as a function of the double root that r has on 
(7r/2,7r). Figure |] shows a plot of r(£) and |m (£)| for the optimal choice of the root. 
Note that the Cohen criterion is clearly satisfied. 

Numerical optimization of the regularity exponent was done by using the Matlab 
Optimization Toolbox. The smoothest wavelet was always found on the boundary of 
the feasible region {v e R 2 : minr v > 0}, that is, the optimal r has a double root on 
(7r/2, 7r). For comparison with the Daubechies family, the most regular scaling function 
with ten coefficients is shown in Figure |5| together with its first derivative and graphs 
of the corresponding Daubechies scaling function. Figure [| shows graphs of the Fourier 
transforms of the scaling functions. 

Table [| lists the locations of the root of m on (tt/2, 7r) for the most regular wavelets. 

The set {v e M 2 : min r v > 0} contains a parametrization of all wavelets with 2N 
coefficients and N — 2 vanishing moments. Together with the results of the previous 
section this shows the regularity exponents in the n z = 1 column of Table |T] are the 
best possible among wavelets with at least N — 2 vanishing moments and 2N > 8 filter 
coefficients. 
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1.8 2 2.2 2.4 2.6 2.8 3 

z 

Figure 3: The Sobolev regularity exponent function of the root z that m and r 

have on the interval (vr/2, it). Here 2iV = 32, n z — 1. 




Figure 4: The function r(£) for 2iV = 32, n 2 = 1 corresponding to optimal regularity. 
The inset shows the graph of |m (£)|. 
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Figure 5: On the top the most regular scaling function and its first derivative, on the 
bottom the Daubechies scaling function and its first derivative. Here 2N = 10. 
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Table 4: Location of the root of r for the most regular wavelets for n z = 1. 
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Figure 6: Graphs of given by fl2|) for the most regular (solid curve) and Daubechies 
wavelets (dashed) with ten coefficients. The blow-up of the interval [4tt, 8tt] in the inset 
shows how the root of mo on (tt/2, tt) improves decay of <f>, hence regularity of 4>. 

4.3 More degrees of freedom 

The approach of the previous section did not work with higher degrees of freedom, i.e., 
with v e M n , n > 3, due to convergence problems with Matlab's constrained optimizer 
(since r is a polynomial of high degree, problems with the constraint min r v > are not 
unexpected). Instead, the optimization was done by considering the Sobolev regularity 
exponent sq as a function of the zeros z = (zi, . . . ,z Hz ) of r on (tt/2, it). The benefits 
are that now we are studying an unconstrained optimization problem and the number of 
degrees of freedom is halved. 

For a fixed value of z we can solve for the parameters v discussed above by requiring 
r v (z) = r' v (z) = 0. This is a linear system, though it may be badly conditioned when r is 
of high degree. (This problem was overcome by using Mathematica and 30 digit precision 
in the computations.) 

The optimal locations for the roots of mo are given in Tables |5|, |6], and |7|. Figure 
shows a typical example of the regularity exponent as a function of the roots and Figure ^ 
shows a plot of r(£) for the optimal choice of the roots. Again the Cohen criterion is 
clearly satisfied. 

As discussed in Section |2] it is still open what happens if r is not required to have as 
many roots on (tt/2, tt) as above. Recall that r is essentially a polynomial (of cosines) of 
high degree. Figure |] shows high peaks between the roots of r. Is it possible to make 
these peaks lower by not requiring r to vanish between them? However, even if this were 
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Table 5: Location of the roots of r for the most regular wavelets for n z = 2. 
the case, it is not clear how the regularity exponent would change. 



5 Other parametrizations 

The parametrization of the wavelets could also be based on the original result of 
Daubechies [§], Proposition 6.1.2], or on a similar result by Wells ||16|| . Both of these 
were tried, but they suffer from the same difficulty as our parametrization in terms of the 
vector v: there is a non-linear constraint that makes numerical optimization for compli- 
cated filters impossible with the software that was used. (These parametrizations could 
be used as the starting point for a parametrization in terms of roots of m , though.) 

Direct parametrizations of the coefficients Ck have been introduced in ||, [TTj, [T2], 17], 
but in contrast to our approach the moment condition ((^) is then a system of non- 
linear equations of the parameters (see in particular flTT| , Theorem 2]). For wavelets with 



2N filter coefficients and M vanishing moments there are N parameters and M — 1 
non-linear equations. The equations can be used to solve for M — 1 of the parameters 
when the rest are kept fixed. Numerical experiments show, however, that typically there 
are multiple solutions that correspond to wavelets with different regularity, making this 
approach unsuitable for optimizing regularity. Moreover, these parametrizations provide 
an expression directly for m (£), but for the computation of the regularity exponent we 
need r(£). 

Only the parametrization in terms of the roots of mo has provided the computational 
simplicity necessary for optimizing regularity when the filter is complicated: recall that 
there are only linear equations to be solved. 
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Table 6: Location of the roots of r for the most regular wavelets for n 
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Table 7: Location of the roots of r for the most regular wavelets for n 
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Figure 8: The function r(£) for 2N = 32, n z = 4 corresponding to optimal regularity. 
The inset shows the graph of |mo(£)|- 
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6 Filter coefficients 

The optimization provides us with an expression for |m (£)|. The (low-pass) filter coeffi- 
cients Cfc in ([I]) are then easy to compute by using spectral factorization p. 172]. The 
coefficients for selected most regular wavelets are listed in Tables ^HlO|. To guarantee that 



mo satisfies the orthonormality condition dJ) to high precision, the factorization was done 
with Mathematica using 30 digit precision and checked by repeating the computations 
with 60 digit precision: the results are accurate to the number of decimals given in the 
tables. 



7 Conclusions 

We have constructed wavelets that are more regular than any of the wavelets that have 
appeared in the literature with the same number of filter coefficients, thus improving 
the ratio of regularity to complexity. Our methods show that a numerical optimization 
approach provides the flexibility necessary to construct highly regular wavelets while 
keeping the number of vanishing moments of the wavelet as high as possible. 

We have also shown our results are optimal for wavelets with 2iV filter coefficients 
and at least N — 2 vanishing moments. 
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Ck 

1.807084186243315e-l 
6.272955371644549e-l 
7.021176047824324e-l 
1.120128480002895e-l 
-2 . 446469866533168e-l 
-2 . 970511286791226e-2 
8.518480911807988e-2 
-7.179751673520435e-3 
-1 . 625706468497944e-2 
4.683260563235789e-3 



Table 8: 2N = 10, n z = 1 



Ck 

4 . 409469394058257e-2 
2.599093138855805e-l 
6.057694016921036e-l 
6.371421499968423e-l 
1.364367462665367e-l 
-2 . 883826797756209e-l 
-1 . 270531514862528e-l 
1.555314597190144e-l 
6.802158361295070e-2 
-8 . 862450100620332e-2 
-2 . 399047021783915e-2 
4.548910378211001e-2 
2.025062066984338e-3 
-1.790050908684461e-2 
3.315675665790194e-3 
4.335970768455452e-3 
-1 . 852310988320677e-3 
-3 . 359209223085526e-4 
3.395506340120817e-4 
-5 . 760617447778257e-5 



Table 9: 2N = 20, n z = 2 



Wavelets with optimal Sobolev regularity 



9.286962261105670e-3 
8.065009178408755e-2 
2.975638977797063e-l 
5.849583107003815e-l 
5.940950288790657e-l 
1.464131039850201e-l 
-2 . 899897992799635e-l 
-1.888671876256031e-l 
1.529318071875427e-l 
1.393659002859078e-l 
-9 . 487951734281382e-2 
-8.688845136891955e-2 
6.416864747737246e-2 
4.586586496411027e-2 
-4 . 242829324588881e-2 
-1 . 870129886474555e-2 
2.492219541071851e-2 
4.323699937678751e-3 
-1 . 195215551785778e-2 
7.996165686836154e-4 
4.259303332169016e-3 
-1.263090636020650e-3 
-9 . 623632092949613e-4 
5.655466161125865e-4 
7.527472879191102e-5 
-1 . 260995097154413e-4 
2 . 020452355886283e-5 
1.026632536894423e-5 
-4 . 41 1797664714925e-6 
5.080242007105260e-7 



Table 10: 2N = 30, n z = 3 



Wavelets with optimal Sobolev regularity 
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